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TIDAL COMPUTATIONS IN SHALLOW WATER 


J.J. Dronkers* 
J.C. Schonfeld* 


SYNOPSIS 


The computation of tidal elevations and currents in shallow coastal waters 
may serve various practical purposes. Mathematically the problem is so in- 
volved that no simple procedure to be generally accepted exists. There are 
mainly three ways of approach which are illustrated in this paper by expound- 
ing the computation methods developed in Europe, in particular in the Nether- 
lands. It depends on the sort of problem to be solved which method should be 
adopted as the most appropriate. 


INTRODUCTION 


1,1. Purposes of computation. 


Hydraulic engineering in maritime waters is confronted with the tidal mo- 
tion. It depends on the extent of the structures to be planned, how deeply the 


engineer will be interested in the tides. 

When he is concerned with local structures of relatively small extent, so 
that there will be no serious interference with the tidal motion as a whole, it 
will be sufficient to collect observational data on the tide as it exists. The in- 
fluence of the structure on the local pattern can then usually best be investi- 
gated by a model of the local situation. 

If however a substantial interference with the movement of the tides is con- 
templated, it will be necessary to investigate thoroughly the mechanism of the 
tidal motion, in order to predict correctly how the intended interference will 
work out. For this purpose the engineer can have recourse to computations 
and to research on a model ofthe whole estuary. 


Technical projects which may entail such a thorough ‘nvestigation are e.g.: 

1. Land reclamation in an estuary. 

2. Safeguarding low countries along an estuary from flooding by storm 
surges. i 

3. Improvement of draining of low countries along an estuary. 

4. Preventing or impeding the intrusion of salt water through an estuary. 

5. Preventing the attack of a tidal current on a shore. 

6. Improvement of a shipping channel in an estuary. 

7. The construction of a shipping canal in open connection with the sea. 

8. The utilization of the energy of the tides. 


Usually a technical project will cover more than one of the above purposes. 


1. Chief Mathematician, Central Research Division, Netherlands 
Rijkswaterstaat. 
2. Chief Engineer, Central Research Division, Netherlands Rijkswaterstaat. 
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The computations provide informations on water levels (of interest for 
height of seawalls, draining sluices, draught of ships), currents (shore pro- 
tection, silting up or deepening of channels, navigation) and energy (power 
plant). They give indications for the execution of works, in particular for 
closing programs of stream gaps. Sometimes computations may have in- 
fluence on the qgperet design of a project, like in the base of the enclosure of 
the Zuider zee 2). where the location of the main dam was altered according 
to insight gained by computations. None of the many works in the Dutch tidal 
waters executed since, was undertaken without the support of tidal computa- 
tions. 


1,2. Nature of the problem. 


Tidal hydraulics in shaliow water deals with the mechanism of the tidal 
motion in estuaries, inlets, tidal rivers, open canals, lagoons and other coast- 
al waters. For brevity we shall hereafter often speak of ‘‘ estuaries,’’ when 
we mean those shallow coastal waters in general. 

The astronomical tides are generated substantially in the vast oceans and 
thence penetrate directly or through border seas into the coastal waters just 
mentioned. The tidal motion in this final stage may be characterized by the 
following properties: 

i. On the whole the tides belong to the kind of wave phenomena called “ long 
waves,” i.e. waves in which the vertical velocities and accelerations are neg- 
ligible. Only the tidal phenomenon known as “ bore’”’ forms an exception (5, 4). 

2. The estuaries are usually so shallow that various effects which are al- 
most imperceptible in deep water such as bed friction and nonlinear distor- 
tion, become appreciable or even predominant. 

3. Shores and shoals substantially impose the direction of the flow of water. 
The estuary may therefore by considered theoretically as a channel or as a 
network of channels. 

4. The inlet in which the tides penetrate is usually so much more narrow 
and shallow than the sea or ocean whence the tides come, that the reaction by 
the inlet on the sea or ocean is negligible or at most a secondary effect. 
Hence the tidal motion at the offing of the inlet may be considered as the given 
source of the motion in the estuary. 


The tides are seldom of purely astronomical origin. They are in particular 
often affected by meteorological conditions (storm), sometimes in a consid- 
erable degree. Since in shallow water the nonlinear effects are strong, the 
deviations from the astronomical tide (the storm surge) can not be well con- 
sidered separately. On the contrary the composite motion resulting from the 
combined astronomical and meteorological forces must be treated as one in- 
tegral phenomenon. Such tides affected by storms will hereafter be called 
‘‘storm tides.’’ 

Tidal hydraulics deals with storm tides as well as with the normal undis- 
turbed tides. 

The tidal motion in a channel can be described by two differential equations, 
the one expressing the conservation of mass (equation of continuity), and the 
other expressing the equilibrium of forces and momentum in the length direc- 
tion of the channel (dynamical equation). The ways of dealing mathematically 
with these equations can be grouped as follows: 


1, Harmonic methods. The composite tidal motion is resolved into har- 
monic components by Fourier series, and these harmonic components are 
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treated separately while terms for their mutual interaction are introduced. 

2. Direct Methods. The equations are subjected immediately to some 
process of numerical integration, e.g. by an iterative process, power series 
expansions, or by converting the differential equations into equations of finite 
differences. 

3. Characteristic methods. The propagation of the tidal waves is analysed 
on the basis of the theory of the characteristic elements of the differential 
equations. 


1,3. Historical survey. 


At the beginning of the development of tidal hydraulics we meet the work in 
England by Airy(l), which dates from the first half of the 19th century. Airy 
treated the tides as periodic waves which he resolved into harmonic compo- 
nents. He demonstrated that, by the nonlinear character of the propagation, 
an originally purely sinusoidal wave is distorted in such a way that higher har- 
monic components are being introduced. 

After the middle of the 19th century de SaintVenant(2) in France approached 
the propagation of tidal and similar long waves from another side. Although 
the theory of the characteristics is not explicitly mentioned, it yet forms the 
mathematical background of de SaintVenant's work. A contribution in this 
field was given likewise by MacCowan(3) in England. 

Full emphasis on the value of the characteristics for defining the propaga- 
tion of tides is laid by the Belgian Massau(5). His work, which dates from 
1900, has attracted less attention from tidal hydraulicians than it deserved. 

Its merits have only been fully understood about half a century afterwards. 

In the 20th century the question of practically computing the tidal move- 
ment in an estuary comes to the fore. De Vries Broekman(7) (Netherlands) 
was the first to point out the possibility of such a computation by a direct 
method of finite differences, and Reineke(9) likewise developed a direct method 
and applied it to German rivers. 

The art of tidal computations received great stimulus by the decision to 
partially enclose the large estuary of the Zuiderzee in the Netherlands. A 
state committee under the presidence of the great physicist Lorentz was en- 
trusted with the investigation of the tidal problems of the Zuiderzee(12), The 
committee followed two ways of approach. 

Firstly Lorentz contrived by an ingenious artifice to linearize the quad- 
ratic resistance in such a way, that the fundamental harmonic component is 
rendered with great accuracy. On this basis a computation method was de- 
veloped to determine the M2 component of the normal tidal movement in the 
channel network of the Zuiderzee(11). The method was used to predict the 
modifications in the tides after the enclosure. 

Secondly a direct method by power series expansion was developed by 
which some computations of storm tides were performed. 

The work of the Lorentz committee proved to be a fertile ground for the 
further development. The quadratic character of the frictional resistance 
encountered by a tidal flow had always been one of the main practical diffi- 
culties for computations. Although Lévy (France) at the end of the 18th cen- 
tury already had put forward the principle of linearization in computing the 
tide penetrating up a river(4) and Parsons (U.S.A.) had ge a treatment by 
linearized equations in his study of the Cape Cod canal r. the real clue has 


3. A more recent American publication is Pillsbury's ‘‘ Tidal Hydraulics”’ 
(1938), which we must leave out of the discussion to our regret, as we have 
not been able to lay hands on a copy. 
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been the principle of Lorentz. The extension of this principle to rivers with 
a fluvial discharge was taken up by Mazure who developed a method to com- 
pute the M2 tidal component on the Dutch rivers(17), 

The next step in the Netherlands was the analysis of other harmonics as 
done by Airy but with the frictional resistance taken into consideration; 
Dronkers(21), stroband(20) and Schénfeld(28) have each contributed to the solu- 
tion of this problem. 

The work of Van Veen(16) may likewise be mentioned in this context, al- 
though it bears not so much on computation methods as on the technique of the 
electric analogue of a tidal system. 

The direct method by power series of the committee Lorentz was made fit 
for tidal rivers by Dronkers(14), In a later stage a4) power series were con- 


verted into expansions by an iterative process - Many tidal problems 
have been analysed more or less intensively by these methods in the course 
of years(27,32), 


In the post-war period Holsters(19) (Belgium) re-discovered the work of 
his compatriot Massau. The method of cross-differences which he developed 
and presented by the name ‘‘method of lines of influence’’ as an approxi- 
mate characteristic method, should in fact be classified as a direct 
method(33) (cf 4,1). 

The method presented by Lamoen(25) (Belgium) is an approximate charac- 
teristic method in which the nonlinear features of the propagation are neglect- 
ed, but in which the frictional resistance is computed correctly. 

A more refined application of the theory of the characteristics was given 
by Schénfeld(28) who demonstrated the value of the characteristic analysis for 
the fundamental discussion of the propagation of the tides. 

The paper deals with its subject as follows: 

First the mathematical formulation of a tidal problem is discussed (Ch.2). 

Next the groups of methods of computation are expounded in chronological 
order (Chr. 3,4,5). In each chapter the most simple method of the group is 
treated in order to demonstrate the principle. Then the more refined methods 
follow. 

Finally a comparative discussion of the methods of computation is given 
(Ch.6). The fields of their application in European, and more particularly 
Dutch practice, are indicated. Moreover a comparison with model research 
is made. 


List of basic symbols 


A cross-sectional area of streambed 
a(ag) depth below water surface ag = A/bg 
B storing area of a section 
b storing width of water surface b = B/J 
bs surface width of streambed bs = @A/Oh 
Cc Chézy coefficient of flow 
c(Cg, c*, velocity of propagation 
F.G characteristic wave components 
g gradient of gravity 
H total head above datum H = h +v2/2g 
h water level above datum 
i(ip, ig, ig) inclination (= slope) 
imaginary unit j2 = -l 
k conveyance of cross-section of e 
streambed Yas 
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length of section 

inertance of section 

inertance per unit length 

discharge (ebb positive) 

discharge per unit width (length) 
linearized resistance of section 
linearized resistance per unit length 
time 

kinetic factor 

velocity of flow 

quadratic resistance of section 
quadratic resistance per unit length 
coordinate along channel (positive in 
seaward sense) 

characteristic wave impediment (wave 
admission 

density 

time of propagation of section T= fk 
relative phase angle of n-th h 

angular frequency of fundamental tide 
complex amplitude of vertical (hori- 
zontal) harmonic tide 

complex tidal impedance (admittance) 
parallel admittance, of section (per unit 
length 

series impedance, of section (per unit 
length 

complex propagation exponent, of section 
per unit length) 

complex n-th harmonic Fourier coefficient 
of H(Q) 


THE BASIS OF TIDAL COMPUTATIONS 


2,1. Schematization of an estuary. 


Most tidal waters have an irregular shape as well in plan as in longitudinal 
or transversal section. Every irregularity like a shoal, isle, groyne, etc., 
has its influence on the local pattern of flow. It would be a considerable com- 
plication of computations if all these local situations had to be considered in 
detail. Fortunately this is generally not necessary since it is possible to com- 
pute a tidal motion accurately by means of a rather severely schematized 
mathematical model, provided this model represents correctly some particu- 
lar condensed characteristics of the estuary. This must be checked if possible 
by analysing well observed tides. 

We confine ourselves here to the case of an estuary or other tidal water 
with such a small width compared to the wave length of the tide (cf 3,1), that 
the tidal flow is directed mainly in the length of the estuary. 

The bed of the estuary fulfills two hydraulic functions. Firstly the bed 
conveys the flow of water in the length direction. Secondly the bed stores 
quantities of water as the tide rises and returns them during the falling tide. 
Not all the parts of the bed necessarily partake to the same degree in the two 
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functions, Parts of the bed (the channels) partake in both functions. Other 
parts however (like shoals, compartments between groynes, dead branches, 
flooded areas, harbour areas) contribute appreciably to the storing function 
but not or relatively little to the conveying function. 

The estuary is now schematized as a channel that conveys and stores, the 
streambed, and adjacent to it regions that store but do not convey. The veloc- 
ity distribution in the streambed is assumed to be uniform. 

The boundary between the streambed and the adjacent storing regions is 
sometimes well defined by the actual situation, e.g. in a river with a dead 
branch. In other circumstances, when there is in fact a gradual transition so 
that the boundary is ficititious, the schematization is nonetheless justified, 
provided the dimensions of the streambed are defined appropriately. Although 
this can be rationalized, it always remains a matter for a good deal of 
experience. 

The conveying cross-section varies with the water level, not only because 
the depth varies but also because there may be parts of the bed such as shoals, 
that are contributing to the conveying function when the level is high, but not 
when it is low. 

If the cross-section of the streambed is further schematized by a rectangle, 
it may therefore be necessary to apply different schematizations for high and 
low levels. This may entail that a storm tide is computed with another 
schematization than an ordinary tide (cf fig. 1). 

A second schematization is necessary in view of the variation of the cross- 
sections along the conduit. For that reason the conduit is divided into sec- 
tions of not too great length. In each section an average cross-section of the 
streambed is defined and the streambed in the section is treated as a pris- 
matic channel with that average cross-section. 

The total storing area B of the water surface in the section, which is a 
function of the water level h, is divided by the length of the section and this 
quotient is considered as the storing width b of the section. This is generally 
different from the surface width bg of the streambed (b> bs). 

Most tidal currents encounter appreciable losses of head by dissipative 
forces. As a rule the losses by friction along the bottom are predominant but 
there may also be appreciable losses by curvature of the channel, by its widen- 
ing and narrowing and by obstacles like groynes, bridge piers,etc. Provided 
the sections are not too long, it is tolerable to merge all losses of a section 
into an equivalent frictional loss of head distributed uniformly along the 
section. 

When an estuary is formed by a network of channels, the shoals between 
the channels may be divided into separate storing regions of the adjacent 
channels, if necessary with a correction for transmission of water over the 
shoals from one channel to another. 

A complete analysis of the estuary must deal with the channel network in 
all its details. For more restricted purposes however, parallel channels may 
be schematized by replacing them by one channel with a composed cross- 
section. 

The admissable length of the sections depends partly on the regularity of 
the estuary and partly on the character of the tidal motion. Ina very regular 
conduit the sections may be longer than in a very irregular one. Even ina 
very irregular channel, however, it is sufficient that the length of the sections 
be small compared to the wave length of the tide. In Dutch tidal practice sec- 
tions of 5 to 10 km are used as a rule. 
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2,2. The differential equations. 


The tidal motion in the length direction of an estuary is mathematically 
described by two differential equations. They can be derived by considering 
mass and momentum. 


Continuity equation. This follows from the law of conservation of mass: 


Oh 
(201) ox + b 0. 
It expresses that the difference in discharge between cross-sections x and 
x + dx,and the accumulation or evacuation of water by rising or falling of 
the level, must balance each other. 
When a supplementary discharge per unit length qj, e.g. to an area that 
is being flooded over a dyke, must be accounted for, we have 


where the term q; may depend on the head of water. 


mical equation. This is based on Newton's law, which is not easily 
applied directly since the mass of water we consider is variable by the trans- 
port between the streambed and the adjacent storing regions. For this reason 
it is more convenient to use the law of conservation and variation of momentum 
per unit length and in the longitudinal direction of the estuary. From this 
equation we subtract Pv times (201) and after introduction of some approxi- 
mations of minor importance, we arrive at the dynamical equation 


203 — +-—- (l- = 
(203) dx’ cat (1-%v + ip = 0 


The first term represents the gradient of the total head and the second term 
the acceleration of the velocity field. The third term accounts for the convec- 
tion of momentum by the water transported to or from the adjacent storing 
regions. When ¥ = 0, there is no convection. This means that for instance 
water emerging from the storing regions and joining the current in the stream- 
bed, derives its momentum entirely from the motion in the streambed. 

The last term of (203), the resistance slope i,, represents the dissipative 
forces which are all quadratic in v so that we may put 
(204) ip = = wIQiQ. 

K 

Here K is Bakhmeteff's conveyance and w represents the quadratic resist- 
ance per unit length, which depends on the water level h. It can usually be 
treated sufficiently accurately as a frictional resistance by using Chézy's 
formula, which yields 


(205) 


The coefficient C is an empirical quantity. 
It must be observed that no Coriolis or centrifugal forces occur in (203) be- 
cause they are irrelevant as far as the flow in the length direction only has to 
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be considered. In computing cross currents over shoals between adjacent 

channels, it may be necessary to account for Coriolis or centrifugal forces. 
In case of a storm tide it may be necessary to introduce a term for the 

forces exerted by the atmosphere, so that we extend (203) as follows: 


(206) Vv +ig +i, = 0. 


Here i,, the wind slope, represents the force exerted by the wind on the sur- 
face and ig represents the barometric gradient. 


Other forms of the equations. A drawback of the equations (201) and (203) 
is that there appear four dependent variables, H and Q as well as h and v. 
Now it is easier to eliminate H and Q than h and v, but unfortunately this is 
of little use since at the transitions between the sections the quantities H and 
Q are to be treated as continuous and not h and v in which jumps are to be 
accounted for. For this reason we shall at least eliminate @v/ Ot and Oh/dt 
by using the relations H = h + v2/2gand Q = Av. Putting y = 0 we 
obtain 


eS, .§ Bese 
Ox 1-v*/y,2 1-v4/y.? at 


(207) 


OH 1+ Bv2/ve2 bmv OH 


Ox 1 - v4/y,2 ot 1-v2/y,2 dt 


Here m = 1/gA denotes the inertance per unit length (cf (35)) and ve = 
gA/bg is the critical velocity (cf 5,3.). Moreover § is put for (b - bg) /bg. 
When v is negligible with respect to v,, the third terms in (207) and (208) 

are small compared to the other terms. Then the following simplification is 

justified: 


@Q OH 
(209) —+b—=0 
ox 


Ox ot 


In this case we may as a rule put m constant and this will often be justified 


likewise with b and w. 
If v/v, is so great that the third terms in (207) and (208) may not be dropped 


still in most cases v2¢¢ v,2, so that we arrive at 


OF 22 - 0 


Ox ot Ot 


(211) 


ry 2UbQ 0; 
where U = 1/2gA2 = 1/2gm? is the kinetic factor (UQ?2 is the velocity head) 
In (211) and (212) we may treat U as a constant as a rule, and b, m and w as 
functions of H, hence neglecting v2/2g in the determination oi these 
coefficients. 
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As the tidal motion is often largely subcritical (cf 5,3.), the equations (211) 
and (212) are sufficiently correct usually. 


Energy equation. In case of designing a tidal power plant the energy equa- 
tion becomes relevant: 


py? (b - bs) 


h 
(pgHQ) + (1/2 pAv2) + pgbh + - 


(213) pgHq; + pg Qj se pgQis + ig) 


This equation can either be deduced directly from the law of conservation and 
dissipation of energy, or by adding pgH times (202) to pgQ times (206). 


2,3. Particularizing conditions, 


When the schematization of an estuary has been fixed and the coefficients of 
Chezy have been determined, a tidal motion in the estuary can be defined by a 
set of particularizing conditions, usually involving boundary conditions (the 
number of which depends on the complexity of the estuary system) and two 
initial conditions in the whole system or the equivalent of them. 

Where an estuary debouches in sea, the tidal motion in the sea generates 
the motion in the estuary. After careful consideration of the interaction of 
the two bodies of water, it is as a rule possible to set up a boundary condition 
for the estuary involving the total head at the offing as given function of time. 

At the landward end of an estuary we have generally a condition involving 
the discharge. At a closed end the discharge is obviously zero,and up a tidal 
river the discharge must approach the fluvial discharge asymptotically. 

The computation of the tidal motion in an estuary must moreover observe 
boundary conditions at every transition between different channels. 

When a channel is continued by another channel of different cross-section, 
it should as a rule be assumed that the total head at the junction is the same 
in the extremities of both channels. The discharge is likewise the same. 
These are likewise the conditions to be imposed at the transition between two 
sections of an estuary where no particular interference with head of discharge 
prevails. If there is a narrow pass or another obstacle between the two chan- 
nels, a loss of head must be accounted for, and when water is discharged to 
or from the junction from aside, e.g. by a sluice, a difference in discharge in 
the two channels is introduced. 

At a junction of three or more channels, the sum of the discharges through 
the channels to the junction is zero. There are moreover conditions for the 
differences in head between the channels. 

Generally there are in total as many conditions at a junction as there are 
channels meeting there, and hence there is one boundary condition per ex- 
tremity of a channel. 

When the heads and discharges at a definite instant are given throughout the 
whole estuary, we can use this as a double initial condition. 

Often it is very difficult to obtain sufficient direct observational data to 
construct such a double initial condition. It is therefore of great practical 
value that other more suitable conditions equivalent to the initial conditions 
are possible, 

Firstly we may consider a purely periodic tide. Then the condition that all 
heads and discharges are periodic functions of the time with a given period, 
replaces the initial conditions. 
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Secondly we may use the fact that the influence of an initial condition on 
the subsequent motion decays and dies out gradually. It is therefore possible 
to compute correctly the tidal motion in an interval of time in which we are 
interested, by starting from inaccurate initial conditions, provided these con- 
ditions lie sufficiently far in the past. The time of decay to be observed de- 
pends on the degree of inaccuracy of the initial conditions and on the proper- 
ties of the estuary system, in particular its extent. 

The boundary and initial or periodicity conditions define the particular mo- 
tion under consideration, which may belong to one of the following types: 


. An ordinary tide on a particular day. 

. An average tide,usually a lunar mean tide, either diurnal or semi- 
diurnal. 

. An average spring tide or an average neap tide. 

. A particular observed storm tide. 

. A hypothetical storm tide, generally of excessive height. 

. Atide on a river encountering a fluvial flood. 


When there are more observational data on a tidal motion than needed to 
supply the necessary particularizing conditions, the redundant data may be 
used for checking. For the Chézy coefficient is so much liable to variations 
and moreover related so closely to the manner of schematizing, that a check 
as mentioned is practically indispensable in most cases. 

In principle the value of C should be determined for each section separate- 
ly and as a function of time. Judging from the computations of Faure(34) for 
the Gironde estuary, the variations in C may then be very considerable. Ac- 
cording to the Dutch practice however, variations in C can as a rule be made 
relatively small by careful schematization, although rather great deviations 
near slack water cannot always be eliminated. This is of little practical con- 
sequence since the resistance near slackwater is weak and hence relatively 
great errors in C are permissable then. For this reason good results are 
obtained by taking C constant throughout large parts of the estuary system 
and throughout the entire tidal period or the entire flood or ebb interval. The 
value of C then varies from about 50 m’ /sec in the shallower rivers to 
70 m’ ? /sec in the deep inlets. 

When a new canal is dug or when an estuary or part of it is modified radi- 
cally, the value of C has to be assumed. It is then recommendable to esti- 
mate the possible deviation of the assumed value and to compute the influence 
of such a deviation. 


INTEGRATION BY HARMONIC COMPONENTS 


In this chapter we confine ourselves to the periodic tide. First the sim- 
plest method to deal with such a tide is expounded: the equations are linear- 
ized which makes it possible to consider the tide as sinusoidal (3,1). Next the 
nonlinear terms are treated and the interaction of harmonic components is in- 
vestigated. The formulae dealing with a second harmonic are developed more 
in detail (3,2 - 3,4). 


3,1. Single-harmonic method. 


Suppose that b, m, and w in (209) and (210) may approximately be put con- 
stant. Then all terms in these equations are linear, except the resistance 
term which is nonlinear in Q. 
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Now consider an estuary without or with little fluvial discharge, where the 
tidal currents vary approximately by a sinusoidal trend. According to 
Lorentz(12) we may then replace the quadratic resistance by « linear resist- 
ance rQ where 


(301) r= =0.85wiQi 
7 


Here | Q | denotes the amplitude of the tidal flow. 

The relation (301) was set up by Lorentz on the assumption that the dis- 
sipation by the ficititious linear resistance should equal that by the real quad- 
ratic resistance. Afterwards Mazure(17) showed that (301) can be obtained as 
well by a harmonic analysis. This analysis can also be applied if there is an 
appreciabie fluvial discharge Qo, in combination with a tidal flow Q; = re Q 
exp jwt. Then we find 


w/Qo + Qy + = + + higher harm., 
where - 2 1 2 
a) = Q| ~ w(1.27/Q| + 0.23 b) ry = + 
and 


(302) a) r = 2wk,|Q| ~w(0.85/Q| + b) = 2wQ, 
when a) Qo< orb) IQI< Q. For 


a discussion of the analysis of 3,3 and in particular (330). 

We can separate the mean motion Q, and the tide Q, (cf 3,4), and here we 
shall confine ourselves to the tide. Then we determine the linear resistance r 
by (301) or (302), using estimations for Q and if necessary Q» to be checked 
afterwards, and obtain the linearized equations 


OH 


(303) Ox 


Ox at 


in which b, m and r are now to be considered as given functions of x. We 
shall for the moment confine ourselves to the case that b, m and Yr are con- 
stants, at least section-wise. In an appendix we shall deal briefly with the 
variability of b, m and r. 

The equations (303) and (304) admit periodic solutions of the sinusoidal 
form 


(305) H= re He™ 1H! cos + arg H) 


(306) Q = re Qe 1Q! cos Wt + arg Q). 
Here H and Q, satisfying the ordinary differential equations 


(307) + jwbH = 0 


(308) + Gum +1) Q=0, 


denote the complex amplitudes of the vertical and horizontal tide, i.e. the 
modulus represents the amplitude and the argument represents the phase of 
the tide. 

Both for the physical discussion and for the practical solution of the above 
equations it is convenient to introduce the tidal impedance Z = H/Q and the 
tidal admittance Y = 1/Z = Q/H (cf(28) ch. Tsect. 23). | ¥! represents the 
quotient of the amplitudes of horizontal and vertical tide whereas arg Y 
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corresponds to the angle of phase lead of the horizontal with respect to the 
vertical tide. From (307) and (308) it can be deduced that Z or Y must satis- 
fy the differential equation of Riccati 

dZ_. 2 aY y? 
(309) a) - or b) = -Yp » 


where y, = jwb and z, = jwm +r. 
The general solution of (309a) or b) is 
Zo - tanh kx Yo - Ye tanh kx 


310 a) Z=———>——__ orb) Y= 


where Z or Y is an integration parameter, whereas furthermore 
k =VYp Zs ; Ze z,/y Ye = = VYp/25 


From any solution Z(x) or Y (x) we can derive solutions for Q and H by 
(307) or (308): 


(311) a) Q= Qo exp-KQ(x) and b) H = QoZ (x) exp-Ko (x) 


or 
(312) a) H = Hp exp-Ky(x) and b) Q= Hy Y(x) exp -K p(x) 


where 


x x 
Ko Z or Ky = Yz,dx 


and where Hp or is an integration parameter. 
From (310a) and (311) or from (310b) and (312) we deduce the general solu- 
tion for H and Q: 


(313) H = Hp cosh kx - Ze sinh kx 
(314) Q = Qp coshkx - YgHp sinhkx . 


This may also be obtained by more conventional methods from (307) and (308) 
or by using the particular solutions to the discussion of which we are now 


proceeding: 
Now we put in particular Yp = + Y, in (310b). This yields the elementary 
solutions Y= Ye and Y = -Y¢ which are constant. From this we derive the 


solutions 
H = Hp exp-kx and Q = HoYe exp-kx , 


for H and Q, from which follows 


H = | Ho! cos (wt - (imk)x + arg Ho) 


-(rek)x 


Q= | Holl Yele cos (wt - (imk)x + arg Hg + arg Ye) 


This represents a harmonic wave with the wave length (27/im k) and travelling 
with the phase velocity (w/im k) in the positive sense of x. The wave is purely 
periodic in t and damped periodic in x. The damping is exponential at the rate 
(re k) per unit length. Hence k is called the complex propagation exponent per 
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unit length. The horizontal tide Q leads by the phase angle (arg Y) with 
respect to the vertical tide H. The solutions 


H = Ho exp kx and Q = -Ho Ye exp kx 
derived from Y = -Y,, represent waves travelling in the negative sense (cf (28) 
Ch. 4, sect. 23). 


The interference of two waves travelling in opposite senses is represented 
by superposition of the corresponding solutions. In this way we arrive at 


(315) H = H* exp-kx + H_ exp. kx. 


(316) Q = Y,H" exp-kx - exp kx. 


Here H+ and H™ are integration parameters. The reader may verify that (315) 
and (316) are an alternative form of the general solution (313) and (314), by 
putting Hp = H+ + H~ and ZeQo = Ht - H-. 

When we consider another particular solution Z(x) or Y(x), we arrive at 
other types of solution for H and Q. By putting Yp = 0 in (310b) for instance, 
and then substituting for Y in (312), we obtain all the solutions for which Q = 
0, at x= 0. Ina similar way Zp = 0 in (310a) yields all the solutions for which 
H = 0 at x= 0. Such solutions may be interpreted as standing harmonic waves 
(ct(28) Ch. 4, sect. 23). 


Appendix to 3,1. 


In order to deal with the variations of b, m and r in dependence on x, we 
divide the estuary in sections so small that in each of them we are allowed to 
take mean values for b, m and r. We may then apply (313) and (314) from 
section to section. This demands much computing labour which often can be 
reduced considerably by making use of the functions Y and Z, in particular 
when we can set up a boundary condition for Y or Z; this is often possible. 
Then we compute Y or Z by (310b) or (310a) from section to section, and 
thence deduce H and Q. 

A slightly different procedure was followed by Dronkers (21), who first com- 
puted the argument of Y by relatively long sections, utilizing the fact that on 
many rivers arg Y varies slowly with x. 

In many cases the sections have to be so short in view of the variability of 
b, m or r, that the integration procedure can be simplified to a finite differ - 
ence calculus. Suppose there is a boundary condition for Z. Then Z is com- 
puted from section to section by finite differences as follows: 

Let Z, and Z,, be the values of Z at the ends of a section (xg, xp) with the 
length § = xp - Xg. Then by (310a) approximately 


(317) Zy - Za - Z, 


where Yp = ypt= jwB is the parallel admittance of the section and Z, = Z,/= 
jwM + R is its series impedance. Moreover Zp, = 1/2 (Zg + Zp). When either 
Za Or Zp is known, Zm can easily be estimated fairly correctly and then a 
construction according to (317) yields Z, or Z, respectively. The estimation 
of Zy, is then checked and if necessary the construction is repeated. 

For numerical computing it is more convenient to modify (317) into 


(318) Zp e Za =YpZ,Zp ry Zs 


from which either Z, or Z, is easily solved when Z, or Z, is known. 
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If Z becomes too great (say |Z! >VZ,/Yp) the variations Zp) - Z, become 
excessive and integration of (310b) is more accurate then. 

In fig. 2a a graphical construction for the Panama sea level canal is repre- 
sented. If the Caribbean Sea were entirely tideless, the boundary condition 
Z = 0 would hold good there. Since there is some tidal motion Ha (which is 
given),Z, is not zero but relatively small. This small value can be computed 
with a very satisfactory absolute accuracy by Z, = H A/Qa; even if we use a 
rather crude estimation for Qa. Such an estimation may be obtained as dis- 
cussed further below. Hence we start from Z, as boundary condition and 
construct Z sectionwise from A to P by (317) and then determine Q by 


Xp 
Q = Qp exp K, where K = z Yp2m 


is computed sectionwise; furthermore Qp= H,/Zp follows from the given 
Pacific vertical tide. Finally H = ZQ in virtue of (311 b). 

Additions or subtractions are performed by vector construction in the 
diagram whereas the multiplications are performed by adding arguments con- 
structively and multiplying moduli by means of a slide rule. 

After having finished the constructions the estimated discharges used in de- 
fining the resistance by (301), and moveover Za, are checked and the compu- 
tation is repeated if necessary. 

The computation was executed for a schematized canal of 72 km length, 

180 m width and a depth below mean level varying from 18 m at the Atlantic 
to 21 m at the Pacific end. Cheay.s coefficient was put 74 m’”? /sec. These 
are the assumptions of Lamoen(25), For comparison the harmonic analysis 
of the results of an exact computation by characteristics of (cf 5,3 and fig. 6) 
are likewise represented in fig. 2. 

A simpie way to estimate fairly correctly the discharges in the canal, is as 
follows: 

Let W, M and B denote the resistance, inertance and storing area of the 
whole canal. Let Qc represent the discharge in the middle C of the canal. 
Now according to (301) we put R = 0.85 W/ Q-/ and then we deduce from (308) 
the approximation 


Ha - Hp= Qe 


where the left-hand member is known. Taking absolute values of both mem- 
bers yields a quadratic equation in | Qc! 2 with a unique solution by virtue of 
| Qclbeing real and positive. After substitution of !Qc!,Qc can be 
solved. 

Then we put Ho = 1/2 (Hy + Hp) and compute with the aid of (307): 


1 1 1 . 1 1 1 
Qa =Q% + 5 Ha Mel. 


These results for Q, which are represented in fig. 2, can be used as basis for 
the above more detailed analysis. 


3,2. Preparations for multiple harmonic methods. 


The approximation of a tide by a simple sine function, however useful for 
exploring a tidal problem roughly, is too crude in many cases when a more 
detailed investigation of the tidal phenomena is demanded. We can then try to 
treat the tidal motion as a purely periodic phenomenon composed of a 
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fundamental and higher harmonic components. The period will as a rule be 
the period of the lunar tide (12 hours 25 minutes). 

The computation of the fundamental component is relatively easy as long 
as the higher harmonic components are not too strong (say less than 40 per 
cent of the fundamental). Then the influence of the latter on the fundamental 
is negligible so that the fundamental may be computed substantially along the 
lines of the preceding section 

The higher harmonic components demand much more computation labour 
and this labour increases disproportionally with the number of harmonic com- 
ponents to be computed, owing to the strong mutual interaction of higher har- 
monic components. This is associated with the fact that the deviations of a 
tidal curve from the simple sine form, or from a combination of a fundamental 
and a second harmonic component, are generally not well represented by one 
single higher harmonic component. 

In regularly shaped rivers the second harmonic component is as a rule a 
fraction of the fundamental, and the third is a fraction of the second harmonic 
(in the Dutch rivers a half or less, and one third respectively). In such cir- 
cumstances a computation of zero, first and second harmonic component only, 
will meet most practical requirements. Occasionally the higher harmonics 
are so small that they may be neglected altogether. 

In other cases, e.g. in more irregularly shaped rivers and estuaries, the 
second and third harmonic are possibly appreciable and of equal order of 
magnitude. It would then be necessary to compute the third harmonic as well 
since the first and third harmonics together produce a second and other har- 
monics owing to the nonlinear terms in the differential equations,in particular 
the quadratic resistance. 

It may be assumed that it is economic as a rule to compute the second har- 
monic. When this is no longer sufficient so that further components are re- 
quired, abandoning the harmonic method for an exact method, e.g. by a char- 
acteristic analysis (cf Ch. 5),is usually preferable. For that reason we shall 
hereafter confine ourselves to developing the formulae for the zero, first and 
second harmonic. The formulae for the third and higher harmonics may be 
derived if necessary along similar lines of thought (cf the appendix to 3,3). 

We base ourselves on (211) and (212) where we treat b, m and w as func- 
tions of H; we put U constant. 

We consider a periodic tide with period @ and fundamental angular fre- 
quency w = 27/@ and expand H and Q in Fourier series. Hence 


H = Ho + Hy + Ho 


where 
H, = H,e"™ +H,e = 2|H,! cos@t + arg H,). 


Here H,, denotes the complex conjugate of H,- The constant H) is the mean 

head, H, the fundamental tide, Hj the second harmonic tide, etc. The modu- 

lus of the complex constant H, represents half the amplitude of the n-th har- 

monic component and arg H, its phase. So H of the preceding section is 2H. 
In the same way we analyse Q: 


Q = Qy + Qy + Qo + Q, +cc (n>1) - 


Here Qo is the mean discharge (on a river identical with the fluvial dis- 
charge). Furthermore cc denotes the complex conjugate of the preceding 
term. 
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Usually higher harmonic components are weaker than lower ones. Yet this 
is not at all a rule without exceptions. We shall assume however, that the 
fundamental dominates over the second and higher harmonic components. 
Then we may usually assume moreover that the variations of b, m and w in 
the course of time are substantially defined by the fundamental vertical tide. 
So we put 


(319) b=b Hy . 


Here 


b = bltig + 2 cos ; 


7 
1 1 
bl) | + cos dé , 


where 


g = wt+arg H; - 


Usually these coefficients b() and b(1) are almost independent of the ampli- 
tude Hy! . 
Similarly we put 


(320) m = m) m(1) Hy . 


where m() and m(1) may be defined in the same way as b(0) and b(1). Instead 
we may apply an analysis as following below for w. 
In a channel with a rectangular cross section we have 


Now neglecting the velocity head we may put 


+H, + Ho 


where a, is the average depth during the tidal period. As Hg is assumed to be 
smali compared to ap + Hy, and H; < 4p, we expand in powers of Hg as follows 


1 3 
YW 
+ H,)° (ag + H,)* 


Hy + 


For the sake of brevity we shall further omit the terms with Hg which are 
often negligible. Then we can put 


(321) w = w(0) wil) elut cc) + (H,2 e2iwt cc) 


where the coefficients w(), w(1) and w(2) are defined by Fourier analysis of 
the factor 


1 
(ay + Hy)? [ag + cos + arg H,)]° 
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This yields 


32,3 
ag? - 4] 


For the general definition of the coefficients c,(¥) cf(28) Ch 14 sec 11. 
It now remains to analyse the quadratic factor |Q/Q inthe resistance. In 
view of the importance of this factor we shall devote a separate section to it. 


3,3. Analysis of the quadratic resistance. 


In case Q keeps the same sign, say +, throughout the entire period, the 
analysis of the quadratic factor |Qi Q= Q? offers no particular difficulties. 
Neglecting third and higher harmonics in Q as well as in Q2, we find by simply 
executing the multiplication of the series for Q with itself: 


(322) Q = (Qo? + 2|Q4|? + 2 2) + [(2Q,Q + elt + cc] 


cc]. 


+ [(Q,? + 2Q99,) 
When Q changes sign during the period, the analysis of 1QiQ becomes 
much more complicated. The first who treated this problem was Mazure(17) 
who confined himself to the case that Q is a simple sine function 


Q = Qo + elt + cc) = Qn +2 cos (wt + arg Q;)- 


Then S = |Q|Q is a non-sinusoidal periodic function which can be decomposed 
in a mean value Sp, a fundamental S, etc. Mazure demonstrated that, in case 
Qo = 0, this fundamental S; is exactly Lorentz’ linearized resistance defined 
by imposing the condition that the linearization should yield the true dissipa- 
tion of energy during an entire period. 

After the work of Mazure it has been tried to extend the theory by consider- 
ing also the higher harmonic components. This encounters great practical 
difficulties however as we explain below: 

In order to perform the integrals of the Fourier analysis the instants at 
which the flow turns have to be determined,because at those instants the factor 
S= |QiQ = + Q2 changes sign likewise. The instants of slack water are de- 
fined by goniometric equations. Even if we neglect the third and higher har- 
monics in Q, this goniometric equation is still equivalent with a quartic alge- 
braical equation and therefore it is not possible to represent the roots by a 
simple formula. Consequently the results of the Fourier analysis can not be 
brought into a workable form either. 

Therefore we must have recourse to approximate procedures. One of these 
procedures consists in approximating the instants cf slack water by the zeros 
of Qo + Q;- This we treat below in connection with Schénfeld's turning 
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w(t) 
C,(3) 
| 


function (cf (28) Ch 14, sect 122). In an appendix we shall deal with more 
refined approximations. 
We introduce a turning function T defined by 


+1 if Q>0 
T(t) = 
-1 if Q<o, 


so that we may put S = |QiQ= TQ. The Fourier coefficients of T, defined by 


7 
-jnut 
T rw dt 
2n 


might easily be computed if we knew the instants of slack water. When there 
are two such instants in a period, we have 

1 
-wty 
(323) a) 


- if ty< t,} 


1 
b) 1, le jnut; iat 


where t, is the instant at which Q turns to the positive and ty} the instant at 
which Q turns to the negative. 
Now we shall approximate ty and t,1 by the zeros of the function Q + Q)- 
We assume Qp< 2/Q,! for otherwise there is usually no slack water at all. 
We suppose Qo> 0 and introduce an auxiliary angle y by 


(324) cos y = Qo/2 1Q,! 


so that 
(325) Q = Qo + Q; = [cos y + cos (ut + arg . 


Hence 
(326) a) wt) =7+y -argQ; ;b) wt)l=7-y-argQ,. 
By substitution of (326) in (323) we deduce 


(327) T = ko" Kn= + ce) , 
where 
Now we proceed further as follows: 


(328) a) =1 b) = (-1)P*2 (n #0). 


(329) TQ? = + Qy + + = + + 2TIQy + Q)) 


+ Qg + + + + 


We confine ourselves, as said before, to the case that Q3 etc. may be ne- 
glected. So we drop the terms 2T(Qp + Q1)Q3 etc. Then the term 
T(Q2 4 .....)2 is likewise negligible as a rule. 
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Applying (326) and (327) and dropping third and higher harmonics yields 
(330) T(Qo + = 4k 1Q4! 2 + | (Qy + cc) + 


+ (Q,2 e2jut + cc) 


(1 
ko = (1 - =) + Cos vy) sin* y 
ky = (1 - 7%) cos y +2 (3 sin y +% sin 37) 


2 1 
kg - +5 sin 2y Sin 47). 


The above result conforms to Mazure's analysis. 
The introduction of the second harmonic of Q produces a number of terms 
of which the following are the most important: 


4ko 


1 
(Q?; Q + ce) 


1-4 -2% 
ko* = (1 7) cos y +> Sin y 


1 2 1 
=5 (1 Sin 2y 


ko 7 Siny g Sin 3y) . 


For the rte of these coefficients we must apply (328) for n = 0, 1, 2, 3, 
and (324) (cf\28) Ch 14, sec 122). 

The above analysis yields fairly accurate results even if the second har- 
monic component is appreciable, say 40 or 50 per cent of the mean and the 
fundamental. This is explained as follows: 

Dropping the second harmonic only affects the approximations for the in- 
stants of slack water. This means that in the interval between the assumed 
and the real instant of slack water, a wrong sign is appended to Q. This is of 
relatively little consequence however, since is small near slack water. 


Appendix to 3,3. 

When the higher harmonics in Q are strong, the above analysis is no longer 
applicable. In this appendix we treat briefly two methods to be considered 
then. 

1. We approximate !Q!Q by a polynomial, e.g. a cubic, as follows: 

Let Qr, + Qy be the greatest and Qr, ~ Qq the smallest value of Q during 
a period in a definite place (Q,>Q,,; otherwise there is no slack water). We 
introduce the parameter p = Q,, /Q, (0< p<1) and put x = (Q - Q,)/24: Then 
we have 


|Q|Q= Q2,|p +x | (p + x). 
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Here 


Now expand | Q |Q in the interval -1<¢ x < 1 by the series: 


IQiQ= S,P,(x) , 


where P,,(x) denotes the polynomials of Legendre. By virtue of the fact that 
these polynomials are normal, we have 


2 2 +1 
=Q (p+x) P, (x) dx. 
a J 


By these integrals the coefficients S, are defined as functions of p. We termi- 
nate after S3 and then obtain the cubic approximation: 
3 


(332) S= 1Q1Q ~ +m + + ng 


The coefficients no, ny, n and n are functions of p (ct?) tig, 105). 

By substituting the Fourier series for Q in (332), the series for S is easily 
deduced. 

An alternative approximation in the form of an odd power polynomial of the 
seventh degree, as deduced by Stroband(29) holds good for the circumstances 
on the Dutch rivers, but the cubic (332) has a considerably wider range of ap- 
plication. 

It has appeared that the procedure by Legendre polynomials is not quite 
free from objections, which make an extension beyond the third degree not 
advisable. For this reason recently the problem has been approached from a 
new angle: 


2. We treat the factor | Q|Q by first analyzing |Q! as follows: 
From the Fourier series for Q we can easily deduce a series for Q2. Then 


we have 


IQi = VQ V Pil +P(t)) , 


n 2n iqut 

P = Q)? 2 and g(t) = 5. (B,e + cc 

Q Q,! nd $(t) qu Pa 
where Bg denotes a set of coefficients depending on the Fourier coefficients of 

If max @(t) during a period is less than 1, we can apply the binomial series 

1 
p=0 \p /‘q=l 4 


In practical applications however, max @may very well be nearly 1 or great- 
er. In that case we write 


(333) 5 Ht) - aft), 


where a is defined as follows: let A be an estimate of max ¢, then we require 


a to satisfy 
vil 
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The value of a is not very sensitive to variations of A. Now by virtue of 


we have 


4 
(334) S= IQiQ = 


by which the Fourier coefficients of S can be deduced. 

The clue of the above method lies in the fact that (333) is most accurate for 
the greater values of |Q/.Perhaps it is less accurate for small values, but 
this is of little consequence since then the product |Q/Q is small. 

In order to demonstrate the value of the above approximations we consider 
an example in which the second harmonic is twice as strong as the fundamen- 
tal: 

Q = cos wt +2 cos 2ut. 
Exact analysis yields 
1Q1Q = 0.36 + 2.62 cos wt + 4. cos 2wt + 0.7 cos 3wt + 0.1 cos 4wt . 
Furthermore we obtain by (332): 

1Q1Q = 0.42 + 2.72 cos wt + 4.17 cos 2ut + 0.9 cos 3wt + 0.4 cos 4ut . 

Finally (334) yields 
1Q |Q = 0.31 + 2.55 cos wt + 4.06 cos 2wt + 0.75 cos 3wt + 0.22 cos 4wt . 


Apparently the latter is the closest approximation. 


3,4. Separation of harmonic components. 


The Fourier expressions derived above are now substituted in the terms 
of the differential equations (211) and (212). Then we have, confining ourselves 
to zero, first and second harmonics: 


0Q 2 wt 
Ox = * dx 


Furthermore (cf (319)): 


at n=l Ot 1 


+ cc) + (2jwb + cc) + (jwb) 2e2jut + cc), 


Here terms which are usually negligibly small have been omitted. In the third 
term of (211) which is small, we confine ourselves to the terms: 


2009 -2Ub() (Qo + 


= -(2jwUb QoQ, et + cc) - (2jwUb Q,? , cc) 


In a similar way we get in the dynamic equations (cf (320), (321), (330) and 
(331)): 
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+ cc) 


-2um'1) | Hy! - 1Q;! sin + Gwm Q, elut , cc) + 


+ (2jwm 9) Qo + cc) - (Gum ") Hy e2jut + cc) 


2UbQ = sin - Qo H, cc) + 


(2jwUb() elt + cc) 


(0) re(Q 2Q.) 


8k, re (HQ; Q,) + 


+ 1 Qy + + 1Q11 Qo + 2 


4k,} wl) | + Cc. 


Here +arg - arg denotes the angle of phase lead of the current 
fundamental Q; with respect to the head fundamental Hj. 

When the above expressions have been substituted for the terms of the dif- 
ferential equations (211) and (212), these equations can be resolved into sepa- 
rate equations for each harmonic component. 

The terms independent of t must satisfy the equations 


(335) dQo 


(336) 


dH 
+ + w(4ub© 2m) sin + 


re ~ 


The underlined terms are as a rule small compared to the main terms which 
are not underlined. Double underlining denotes the smallest terms. 
The coefficients of the factor eJ/t must satisfy 


(337) + - 2jwubM = 0 


(338) 


jum%Q, - + IQ - 2 + 
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+ 4k; Qo + 4k, w'2) - - 1Q41 = 0. 


The underlinings again denote orders of magnitude. 
The coefficients of e2t must satisfy 


(339) dQ 
+ 2jwb) Ho + H,2 - 2jup y Q,2 =0 


(340) <2 — + 25m - ju(m™) + 20 - + 


+ 4kqw()Qy2 + + 2 - 0. 


All these terms are small compared to the main terms in (337) and (338). 

The solution is searched for along the following line: 

First we neglect the double underlined terms in (337) and (338) and solve the 
fundamental tide substantially as described in 3,1. 

Then we drop the double underlined terms in (336) and compute Ho by nu- 
merical integration. Here we can be supposed to know Qp and one boundary 
condition for H (estuary or maritime river), or we have two boundary condi- 
tions for H (canal between two seas or the like). 

Next we substitute the results for the zero and first harmonics in (339) and 
(340). These equations are linear in Hz and Qz and nonhomogenous. They are 
solved by applying the theorem that every solution can be expressed as the 
sum of an arbitrary particular solution and a complementary function being a 
solution of the homogenized subsidiary equations. An arbitrary particular 
solution is easily constructed by integration by finite differences from section 
to section and complementary functions can be determined substantially as 
described in 3,1. 

Finally we correct the fundamental and zero harmonics for the double 
underlined terms. 


We conclude by making two remarks: 

The influence of the small terms presents an intricate question. In order 
to justify the neglect of certain terms, it is not sufficient to verify that each 
of these terms is small. It may be that a rather great number of small terms 
all have the same sign so that they accumulate. If this occurs, it may be worth 
while to compute these terms, at least some of them, in order to get an idea 
of the tendency of their influence. 

If other terms than those presented above have to be introduced, e.g. by 
making use of one of the analyses of the appendix to 3,3, the derivation of the 
formulae follows substantially the same line. 


DIRECT INTEGRATION 
First the finite difference methods are discussed, both in the original quad- 
scheme and in the more recent cross-scheme (4,1). Then the principles of 


the more refined methods of power series and iteration are expounded (4,2). 
This is followed by various applications of these methods (4,3 - 4,5). 
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4,1. Finite difference methods. 


Among the first who proposed the numerical computation of tidal move- 
ments are de Vries Broekman(”) and Reineke'9), Both put forward a direct 
integration of the differential equations by finite differences. More recently 
Holsters(19) has introduced another method which we consider as a modified 
application of the idea of finite difference integration. 

The efficiency of this kind of methods is greatly improved by arranging the 
first order differences symmetrically such that the second order differences 
cancel. In order to discuss this, let the time be divided into relatively small 
intervals of equal duration At whereas the estuary is divided into rather short 
sections. Then a grid is formed in the tx-diagram as illustrated by fig. 3, in 
which AF may represent an estuary with the inlet at F and its landward ex- 
tremity at A. 

Two main procedures are to be distinguished: working by quad-differences 
or by cross-differences. 


Quad-differences. Let Q be expanded in a Taylor series in the environ- 
ment of the centre M of the rectangle 11-21-22-12. Then 


= Qu +5 Ax - Oat At + At? 


2 
12 At 


1 


1 1 2 


2 1 1 


1 


where the index x denotes a differential quotient with respect to x and the in- 
dex t one with respect to t. We deduce approximately 


(401) Qi1 + Q2i - - Q22 | 
2hi2 


Here the first neglected terms are of the second order in Ai and At as com- 
pared to Q, hia . Ina similar way we have 


(402) hy hai + hgg - hy, - h 
2At 
Moreover 


hy, + hyg + hgg + hay 


(403) 
4 


and b = b(h)may be considered as accurate in the same degree as (401) and 
(402). Substitution in (201) yields 


1 
(404) 5 + - + (hay + hog - - = 


where B = bl denotes the storing surface of the section considered. 
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In a Similar manner we reduce (203) to 


1 
(405) Ban - - +; + - Yig) + 8H, = 0, 
where we put 
(406) Hy = Wy2 (h) + 


by Wj2 denoting the resistance of section AB. In (405) the third term of (203) 
has been neglected. 

Similar equations are deduced for the quads 12-22-23-13, 13-23-24-14, etc. 
The whole set of equations is soluble provided a sufficient number of boundary 
and initial values is supplied. 

The method presented is substantially that of deVries Broekman. The meth- 
od of Reineke differs from it in that the drivatives with respect to t are not re- 
duced to differences: they may be determined graphically for instance. 

The method of quad-differences has not been applied very often in practice. 
This is due to the fact that other methods offered greater possibilities. 

Firstly the power series or iterative method (4,2) had advantages with re- 
gard to the accuracy. An improvement of the accuracy can be attained either 
by considering shorter sections, or by introducing higher order corrections. 
The former entails revising the schematization which is very laborious, and 
for the latter the power series and iterative method are better suited. 

Secondly the method of quad-differences is not well suited to deal with 
single boundary conditions at both ends of a channel which forms a problem of 
great practical importance. Such conditions are usually introduced ina 
method of quad-differences by some process of trial and error. On this ac- 
count a method of cross-differences or a characteristic method is often pre- 
ferable. 


Cross-differences. Often we encounter a problem in which either H or Q 
is a given constant or a given function of t at each end of a channel. In those 
cases a method of cross-differences may be appropriate. 

We neglect once more the velocity head v2/ 2g and treat A as a constant and 
B and W as constants or given functions of t. We expand H in the environment 
of point 22. Then: 


1 2 
Ho) = Ho2 + Hy + 2 Hyx 412 + 


1 2 


Suppose that =493 or - is small compared to or Log of the 
same order as Hx, $9 is small compared to Hy. Then 


(407) H ~ Hai - B23 __ 
thos 


is correct, but for terms of second and higher order. In the same way we find 


(408) 


By introducing the approximation | Qj9/Q39 for |Q99 | Qg2 the dynamical 
equation yields 


M2 
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Q32 - Qi2 
f 
|| Qt 


| 


where M2 = Mj2 + Mg3 and W 
In a similar way we apply cobfimuit? uation to point 33: 


(H43 - =0. 


(410) Q =3 
32 934+ 


Here B3 = Bo3 


Now suppose Q(t) is given in Aand H(t) in F (fig 3). We suppose moreover that 
Q12, 214, Haz and Hos are given as initial conditions or may be assumed as such. 
Then we know H91, H93, and Q;9 and compute Qs by (409) In the same way 
we deduce Q34 from Hg3, Hg5 and Q;4. Next we *fnow Q32, Q34, and Hg3 and 
we may compute H4, by (410) and fatvermore H45- 

Then we start anew from the values H4;, H43 and Qg9 in order to compute 
Q52 by applying the dynamical equation to the lozenge 32-41-52-43. Proceed- 
ing in this way we compute alternately sets of Q values by (409) and sets of 
H values by (410). The computation gives the vertical tide in B and D and the 
horizontal tide in C and E. When we want to know the other tides, supplemen- 
tary computations are necessary for instance with the aid of quad-differences. 

Once the boundary and initial conditions given or assumed, the computation 
proceeds without any trial and error. The influence of a boundary condition 
is then introduced step by step into the solution along the sides of the lozenges 
(* lines of influence”), This is decidedly an advantage over a method of quad- 
differences, in which a boundary value at some instant t, must be introduced 
simultaneously in all points of the grid with the abscissa tp. 

The above method, in which the continuity and the dynamical equation are 
applied alternately, can only be used if each boundary condition is either one 
in H or one in Q. This restriction can be removed by applying both the con- 
tinuity and the dynamical equation in every lozenge. It is moreover necessary 
then to apply these equations in the border-triangles 16-25-36, 21-32-41 etc. 
Owing to the lack of symmetry in these triangles, the second order terms can- 
not be made to cancel, so that second order differences must be introduced in 
order to mpintein the standard of accuracy. 

When v 2/2¢ and the variations of A are no longer negligible, their compu- 
tation likewise requires the simultaneous application of both equations in 
every lozenge. 

The method of cross-differences has a certain appearance of affinity to the 
characteristic method (cf Ch. 5), in particular when the section lengths and 
time intervals are chosen such that the “lines of influence” coincide with the 
subcharacteristics. This is not essential at all however, and except when B 
as well as M is constant, it is not even practicable. Let 23-S and 21-S be two 
subcharacteristics. Then the values of H and Q along the segment 21-23 de- 
fine the solution in the triangle 23-S-21. Hence, if point 32 lies within this 
triangle, we may compute the solution in 32 from the data in 21 and 23. When 
32 lies outside the triangle, the computation outreaches the propagation of the 
tidal motion’ In that case phantom waves appear in the solution, mainly with 
the period 4At, which tend to grow, in particular near slack water, because 
then there is but little friction to damp them. The phantom waves remain suf- 
ficiently smal when the lozenges keep within the subcharacteristic triangles 


(cf Holsters( 33)). 
4. Physically the influence of a boundary condition proceeds with the charac- 
teristic velocity of propagation. For that reason we prefer to retain the 
name“ lines of influence” for the subcharacteristics, and accordingly we 

prefer to call Holsters’ method, a method of cross-differences. 
714-28 


4,2. Power series and iterative methods. 


In order to compute storm tides in the Zuiderzee, the Lorentz Commit- 
tee(12) developed an integration method a | power series. This method was 
extended to maritime rivers by Dronkers(14) who later converted it into an 
iterative process 21,24) 

We assume once more a division of intervals of t and sections of x (cf fig. 
3). Now consider the motion in the section AB at the instant t}. The section 
is supposed to be so short that it is admissible to treat b, m, wand U iz (211) 
and (212) as constants, introducing for them the mean values in the section 
AB at the instant tj. 


Power series method. Let the origin for x be chosen at the place A and let 
Ho(t) and Qo(t) be the functions H and Qin A. Then we can expand H and Q in 
the environment of the segment 15-16 (fig. 3) by the MacLaurin series: 


‘ 
(411) H(x,t) = Ho + ( <0 


+ 


2 
(412) Qlx,t) = Q +5 x2 


When Ho and Qo are known, all the other coefficients can be deduced by means 
of the differential equations. Then H({gt;) and Q(f3t;), i.e. H and Q in B at the 
instant t;, are defined provided the series converge. 

The convergence of the expansions in the form (411) and (412) is very dif- 
ficult to be ascertained, both theoretically and practically. This objection has 
been overcome by approaching the problem somewhat differently by an itera- 
tive process. 


Iterative method. We write (211) and (212) in the form 


+ wQ2 + 2bUQ 


then substitute the approximations Q = Qo(t) and H = Ho/(t) in the right hand 
members, and integrate from 0 to x. This yields the first subsequent approxi- 


mations: 
(415) = Qo - bHpx + 2bUQQox 
(416) Hy = Ho - mQox + wQp2x + 2bUQgHox . 


Here the points denote derivatives with respect to t. 
The second subsequent approximations are found by substituting (415) and 
(416) in the right hand members of (413) and (414) and integrating once more 


(417) Qr = + bmQpx2 + bwQyQox2 + Qs (t,x) 


(418) = Hy + bmHox? + bwQoHox? + + (t,x) 


Here Q, and H, denote sets of terms of lower order of magnitude. 
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Continuation of the process yields Qyyj and Hy, etc. Evidently the formulae 
become more and more involved. As a rule the approximations Qy and Hyy 
are sufficient for practical computing. When it is necessary to compute Qy7 
and Hyyy, we must generally as well consider the variations of the coefficients 
b, m etc. 

Generally we have 


x. 
(419) Qn = Qo -b + 2bU Qn-1 Qn-1 ax 
0 0 


0 0 


These functions are polynomials in x. 
It remains to be investigated whether the iterative process is convergent 
or not for x<¢ f, i.e. whether the series 


Q= - , H=Ho +2, (Hy - Hy_y) 


converge (Q= lim Q, ; H= lim H,). 
« n— oo 
Since the terms with U in (413) and (414) are small with respect to the 
other terms as a rule, we shall leave them out of consideration for the mo- 
ment. Then it can be proved by induction: 


3(2k-1) 
= pol, (t)x! 
f=2k-1 f=2k-1 


k+1 


(t)x" - Q b 
2k 2k-1 =2k 


The proof of the convergence follows by a more detailed research of the poly- 
nomials of (421). 

The process is found to be convergent for a group of functions Qp and Ho. 
Owing to the quadratic terms in the differential equations, this group is more 
restricted than if the equations were linear. 

In many tidal problems we may assume that all the derivatives of Qo(t) and 
Ho(t) are finite with respect tot. In this case there is convergence for limit- 
ed values of x, depending on the coefficients of the oom equations. We 
renounce the detailed treatment of this question (cf 4 ). 

The terms derived by iteration are the same as those appearing in the 
power series. By their different grouping however, the convergence can 
more easily be proved, The terms are moreover more easily interpreted 
physically as contributions to the balance of quantities of water or to the bal- 
ance of forces and momentum. This facilitates the quantitative appreciation 
of the terms. 

When Q, and H, at the upper end of a section are known, we can compute 
Qp and Hp at the lower end, or inversely we compute Q, and H, from Q, and 
Hp- 

Hence if we know Q and H at the inlet of an estuary or at the mouth of a 
river, we can compute these quantities up the estuary or river from section 
to section as far as the closed end or until the range is so small that further 
calculation serves no purpose. 
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However in various practical applications we do not know Q(t) at the inlet. 
We can then arrive at this Q by trial and error, checking by means of a bound- 
ary condition at the landward extremity which is usually known. 

Other applications will be discussed in the subsequent sections. 


4,3. Computation of currents from vertical tides; single sections. 


In general it is simpler to record water levels than currents. So the prob- 
lem arises to calculate the discharges and velocities in a section as a function 
of time when the water levels at both ends are known. This can be performed 
by means of the formulae of the preceding section. 

We assume the Chézy coefficient C to be known. Let the section be so 
short (in the Dutch practice usually sections of 10 km or less are considered) 
that the second iteration (418) is sufficient. Putting for x the length of the 
section, Hy; becomes the known vertical tide at one end whereas Ho is the 
known tide at the other end of the section. Thus (418) becomes a nonlinear 
differential equation in Qo of the first order. 

Notwithstanding its nonlinear character it is very simple to determine Qo 
as a function of time numerically, when we know or may assume an initial 
condition for Qo (cf. 2,3). We may do so by converting the derivatives with re- 
spect to t into quotients of finite differences. Then we compute Qp from in- 
terval to interval starting from Qo(t,;) at the instant t,. In practice, intervals 
of about a quarter of an hour usually will do. 

When the vertical tides Hp and H, are periodic it can be shown that the dif- 
ferential equation in Qo has one unique periodic solution. This solution is 
stable and hence, if we start from an arbitrary value Qolty), the integral curve 
will approach the periodic solution asymptotically with increasing time (cf 
fig. 4). 

If we are so furtunate as to know the local vertical tides along the estuary 
at distances of about 10 km or less, we are able to check the schematization 
used. Then we calculate the discharges at the beginning of each section by the 
above method. These discharges however must also satisfy the equation of 
continuity (417). This provides a check of the schematization of the estuary 
and of the Chézy coefficient. 


4,4. Computation of currents from vertical tides; combinations of sections, 


The above method is not well applicable when we only know the vertical 
tides at the ends of a channel with a length largely exceeding 10 km. Then we 
may consider the following method: 

Suppose the channel is so long that it should be divided into two sections. 
In each of these sections we adopt mean values for m, b and w as outlined in 
4,2. We denote by Q;, H; the discharge and head at the beginning of the first 
section, by Q;2, Hyg halfway the section, and so on as fig. 5 shows. 

Now we apply (416) to the centre of a section omitting the secondary terms 
involving U, substitute x = x = fo, x = -1/24;, x= 1/2 and then deduce 


- 2 2 ° 
(422) a) H2 = Hy + Wy Qi 241 -mj Qi2hi b) H3 = Hg w2 Qo3h2 - mg Qo3h2 


Moreover we may put 
1 1 
(423)a) = Q2 +5 = - Haghe - 


Now we put approximately 
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Fig 5 Combination of two sections 


lunerhours 


NG 


40 10". sec 
‘eagle - 
iterative, approximate —++— characteristic, approximate 
iterative, exact — —-— characteristic,exact 


Fig 6 Panama sea level canal Computations by various methods 
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a? 3 
Hg ; Ha3 = Hy +243, 


1 


Hy2 = 5(Hy + Hg) =7 


so that 
1 . 1 


Substituting in (422) we find after some calculation: 


(424) Hg - Hy = + (wy hy + who) Q,? ~ (my 4; + mghg) Qo + 


2 
{ (3Hy + Hg) - bo (Hy + 3H) } . 


Thus we have derived a differential equation for Qg, which we may solve as 
treated in 4,3. After calculating Q; we determine 


(425) = + by 


as a first approximation of Q). 

Starting from Q; and Hj we compute Qo and Hg, and then Q3 and Hg with 
the aid of (417) and (418). Now Hg has to be identical with the known function 
Hg, but there will generally be deviations. Then we can determine a closer 
approximation for Q; by putting for H;9 and Ho the functions computed from 
Q; and Hj. This yields a formula analogous to (424) from which we deter- 
mine Q, again. 

We can apply this method for relatively large channels. In the Dutch 
estuaries the total length may be up to 30 km. 

In order to show an example of the preceding method, we have calculated 
the tidal movement in an open unregulated Panama sea level canal using the 
schematization of Lamoen(‘25) (cf 4,1). We have divided the canal in 4 sections 
and computed a first approximation by the above method. Then refinements 
have been computed by (417) and (418). Both solutions are represented in 


fig. 6. 


4,5. Application to the planning of the enclosure of a tidal river, 


As stated in 1,1, tidal computations form a valuable support in the planning 
of the enclosure of a tidal river. In order to illustrate this, we shall briefly 
discuss the computations with regard to the enclosure of the Brielse Maas, 
which is represented schematically in fig. 7. 

Between A and B the river had to be closed by a dam, and likewise between 
Dand E. C is almost halfway Band D. The sections BC and CD are each 
about 11 km long. The vertical tides at A and E, Hp and Hq are known. The 
distances A-B and D-E are so small that we may put 


Qo = and Q3 = 


As the construction of the dam between A and B progresses, the flow 
passes through a narrowing gap. Between C and D the situation is similar. 
In the gap A-B there is a loss of head 
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Fig 7 Scheme of river Brielsemaas 


Fig 8 Characteristic construction diagrams;a 


Diagrem of itineraries, b Diagram of states of motion. 


Ho - Hy = — Wi Q" » (- for flood) 


mg 1,2 


in case of ebb or flood respectively. Here A, is the cross-sectional area of 
the gap, Ap that area above, and A, that below the gap, whereas n; and Ng are 
coefficients of the gap. 

Likewise 


+ 2 
(427) Hg - Hy W,9, 


holds good for the gap between C and D. 

The basin of the Brielse Maas, represented by the sections BC and CD, 
is treated in the way of 4,4. In (426) and (427) we substitute (425) and an 
analogous expression for Q3. Hence we obtain three equations, (424) (426) 
and (427), by which we may determine Qo, H, and H 

After elimination of Hj and Hg with the aid of (426) and (427) we get a non- 
linear differential equation for Qg which may be solved in a similar manner 


as treated in 4,3 and 4,4. 714-34 
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The following questions concerning the closing of the two gaps were put: 

1. Which gap has to be closed first ? 

2. Is it possible to narrow the gaps in coordination with each other in sucha 
way that the closing of one of the gaps would become easier ? 

3. Which are the values of the velocities in the gaps during the process of 

narrowing ? Both the maximum velocities during a tide, and the slack water 

conditions are important from the engineering point of view. 


The results of the tidal calculations showed that it was preferable to 
narrow first the gap at the seaside. Then the velocities in the other gap 
(riverside) would decrease because the penetration of the tide into the channel 
is obstructed. In fact the gap at the riverside could be closed without dif- 
ficulty. After that the gap at the seaside had to be closed, which was effected 
by sinking a large pontoon at slack tide. 


INTEGRATION ALONG CHARACTERISTICS 


The principle of the characteristics and their bearing on the phenomenon 
of propagation can be most clearly discussed in connection with linear equa- 
tions without resistance (5,1). Next the amendments to be made in order to 
deal with the resistance are expounded (5,2). Furthermore the variability of 
the velocities of propagation is deait with (5,3). Finally shock wave conditions 
are considered (5,4). 


5,1. Elementary theory of propagation . 


We start from (209) and (210) where we consider b and m as constants. In 
this section we make moreover abstraction from the resistance. 


(501) 


(502) 


It may be observed that, according to the disregard of v2/ 2g presumed in 
deriving (209) and (210), H may be interpreted arbitrarily as the water level 
or as the total head. 

The tidal motion is mathematically described by H and Q as functions of t 
and x. This can be represented graphically by using a HQ - diagram in con- 
nection with a tx - diagram (fig. 8). 

Let us consider the point t,, x, in the tx-diagram (A in fig. 8). This means 
that we fix our attention to the place x, at the instant t,. Let H, be the head 
and Q, the discharge in this place at that instant. Then the point (H,, Q,) in 
the HQ -diagram (A' in fig. 8) representing the state of motion in x, at tg, is 
associated with the point (t,, x,) in the tx -diagram (A). 

Any function F of H and Q may be interpreted as a property of the state of 
motion. If a definite value F, of such a function F is observed at the instant 
ta in a place xq, and shortly after at the instant t, in a slighly further place 
Xp, we shall say that the property represented by F is propagated from x, to 
Xp during the interval t, to t). 

The further mathematical development of this idea consists in trying to 
deduce from (501) and (502) one or more equations of the form 


714-35 


oH 
at =0 
OH, 22-5 
Ox ot 


0 
ot Ox 
For the elaboration we may refer to(28) Ch 2. Asa result of the analysis 


there appear to be two such equations which we can find by adding 1/2 Zoco 
times (501) to 1/2 cg or to - 1/2 Cp times (502). The two functions which are 


propagated in the above sense are 


1 1 1 1 
(503) a) F 2 and b) 2 


satisfying the equations 


(504) a) OF oF aG 6G 9 
— --=0 db —=0 
ot Ox ot “0 Ox 

where 

(505) a) cg = 1/¥bm b) Zo = = Vm/b - 


In a geometric point moving in the positive x -sense with the velocity 
dx /dt = - Cy, we have 


OF oF OF OF 
(506) = By dt + = + co) at 


which is zero by virtue of (504a). Hence F retains its value in the moving 
point. This point, which we shall call a wave point, can be said to convey that 
particular value of F. Ina similar way a wave point moving in the negative 
sense with the velocity dx /dt = - cp, conveys a particular value of G. 

The functions F and G may be considered as new coordinates in the HQ- 
plane. Hence a state of motion is as well defined by F and G as by H and Q. 
We shall call F and G the characteristic wave components of the tidal motion. 
The values of these components are conveyed by the wave points. Thus the 
component F is propagated in the positive sense and the component G in the 
negative sense. The velocity of pro tion is cg or -cg. 

The itineraries of the wave points are represented by straight lines in the 
tx-diagram, called subcharacteristics, with slopes cg tol or -cg tol. The 
HQ -points associated with the points of a particular subcharacteristic, all 
correspond to the same value of F or G, and hence are situated on straight 
lines in the HQ-diagram, called contra-subcharacteristics. When the slope of 
the subcharacteristic is cg to1, then by (503a) the slope of the associated contra- 
subcharacteristic is -Zo to 1 (compare sj and cq in fig.8). When the slope of 
the former is -cg to 1, then by (503b) that of the latter is Zg to 1 (sg and co in 
fig. 8). A subcharacteristic and the associated contra-subcharacteristic to- 
gether will be referred to as a characteristic. If the solution H(t,x), Q(t,x) is 
represented by an integral surface in a HQtx-hyperspace, the characteristics 
are curves on that surface and the subcharacteristics and contra-subcharac- 
teristics are projections of the characteristics. 

Let us consider a travelling wave invading a state of rest (fig. 9). The 
water surface in rest having been adopted as zero level, we have H = 0 and 
Q = 0 and hence F = 0 and G = 0 in the undisturbed region. Any wave point P 
moving with the wave and conveying a particular value F,, of F, meets wave 
points coming out of the region of rest and thence conveying the value G = 0. 
So we have from (503a) and b 


(507) a) Hp = Fp and b) Hp = 2 Q, 
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in the wave point P independent of t. This means firstly that there is a defi- 
nite relation (507-b) between H and Q and secondly that the whole configura- 
tion of elevations and depressions and the associated currents, displaces with 
the velocity co of the wave points in virtue of (504-a). 


Fig. 9. Travelling (tidal) wave vertical dimensions greatly exag- 
gerated compared to horizontal dimensions. Wave height exag- 
gerated compared to depth 


Fig 10 a Construction to account for resistance 
b Construction to account for variations of velocities of propagation 


The elevations of the water level are attended by currents in the sense of 
the propagation and the depressions by currents in the opposed sense. The 
points R and U where the original level is attained,also mark the places of 
slack water (cf (507b)). 

It can be observed that in a canal with a rectangular cross-section of width 
b and depth ag, we have m = 1 /gagb and hence 


co = » 


a well-known formula for the velocity of propagation of a travelling wave. 

Apparently the travelling wave character of the water motion is associated 
with the wave points moving with the velocity cg, which we shall call the active 
(manifest, cf(28) Ch. 3, sect.11) wave points. The wave points moving in the 
other sense and conveying the rest value G = 0 are inactive (latent). In the 
state of rest all wave points are inactive. This does not imply however that 
those wave points are at rest themselves. 

Now we consider a wave motion in which the wave points moving in both 
senses are active. In order to define the line of thought, suppose we know the 
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State of motion in a place x; at an instant t; and likewise in xg at tp. So we 
consider the points P and Q in fig. 8a. Let P’ and Q' in the HQ-diagram 
represent the associated states of motion. Now we consider a wave point 
moving with the velocity cp in the positive sense, passing by the place X, at 
the instant t; on its way to the place xg. The itinerary of this wave point is 
the sub-characteristic sj. In the same way Sq is the itinerary of an other 
wave point moving with the velocity -cg in the negative sense and passing by 
the place at the instant 

In the tx-diagram we see that the two wave points meet in a place X3 at the 
instant tz. This event is represented by the intersection S of the subcharac- 
teristics sj and S89. The associated state of motion is easily constructed in 
the HQ-diagram. For the HQ-point S' of this state of motion must lie as well 
on the contra-subcharacteristic cj associated with the first wave point, as on 
the contra-subcharacteristic cg associated with the second wave point. 

Now suppose the head H in the place x4 is controlled by a boundary condi- 
tion. This place for instance is the inlet of the estuary. Then the head Hy, at 
the instant t4, at which the wave point first considered above arrived in the 
place x4, is given. The tx-point R represents the arrival of the wave point. 
The HQ-point R' representing the associated state of motion, must have the 
ordinate Hy. The point R' must moreover lie on the contra-subcharacteristic 
cj. These two conditions define R' entirely. 

If the boundary condition controls Q or a definite function of H and Q, the 
state of motion point R' is defined in a similar way. 

Let now the head in the place x4 be suddenly varied at the very instant of 
the arrival of the wave point so that H becomes H4 + 6H4. The HQ -point rep- 
resenting the state of motion immediately after the variation of head, cannot 
but lie on cy as well as R'. This means that the variation of head 6 H, is at- 
tended by a variation 6Q, = -6H,4/Zo (cf (503a)). 

Any sudden variation of head imposed on the canal (x < x4) provokes a pro- 
portional reaction in the discharge and inversely. The increase of head per 
unit decrease of discharge, Zp, is called the characteristic coefficient of 
wave impediment. The reciprocal Yp is called the wave admission. 

In the same way it is argued that sudden variations imposed in the other 
— (say to a canal with x >x4) provoke likewise reactions such that 6H = 
Zo 6Q. 

o The above constructions illustrate some of the main procedures applied in 
the characteristic approach of tidal problems. Other such procedures serve 
to deal with reflections at widenings or narrowings of a channel, at junctions 
of channels etc. 

Before we set out to describe more systematically how a tidal computation 
by characteristics is performed, we must first pay attention to the resistance 
in view of its practical importance. 


5,2. Influence of frictional resistance on propagation. 


We add again 1/2 Zocp times (209) to 1/2 cg or -1/2 Cp times (210) so that 
instead of (504): 


OF OF 1 0G eG 1 
(508) a) + + = 0; b) ry cgw|Q/Q 


where F and G are still defined by (503). The term w Q|Q may be consid- 
ered as a function of F and G. If w is supposed constant, then by (503) 
w 1Q|Q=wYo?|F -G| (F - G) 
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formulates the resistance as function of F and G. 
The function F in a wave point moving with the velocity cg now no longer 
preserves its value. Substitution from (508-a) in (506) yields 


(509) dF = - cow dt. 


Hence F decreases gradually if Q is positive and increases if Q is negative. 
Similarly G in a wave point moving with the velocity -cg increases if Q is 
positive and decreases if Q is negative. 

By the above mathematical arrangement it is possible to explain the be- 
haviour of tidal and similar waves in terms of propagation and attenuation as- 
sociated therewith. The following may illustrate this: 

Consider once more the wave of fig. 9. In the region of rest we still have 
F =0andG=0. Nowa wave point emerging from that region to meet the 
wave, conveys the value G = 0 until it encounters the foot F of the wave. Then 
it enters a region where Q>0. Hence the value of G will begin to increase and 
G becomes positive. So the encountering wave point is activated. Since gen- 
erally by virtue of (503-b) 


(510) ZpQ=H- 2G, 


the currents will be weaker than if there were no resistance. 

By consequence of (510) the place of slack water will be found at a positive 
elevation H = 2G2 0, say in S. At this place G in the receding wave point 
reaches its maximum value and G will begin to decrease behind S. The 
strongest currents are by virtue of (510) found where 


OH , 
6x Ox 


Since 0G/Ox is negative between F and S, the place searched for must have 
OH/Ox negative as well. Hence the maximum flow no longer coincides with 
the top T, but with a more forward point Q. 

The assumption of a region of rest ahead of the wave is not essential for 
the conclusions, although for other assumptions the line of argument is slight- 
ly more complicated. Hence it is explained why in a sine wave subject to re- 
sistance, the horizontal tide has a phase lead with respect to the vertical tide 
(cf 3,1). 

When the line of thought is followed up further it can be explained that the 
displacement of the top of a travelling wave.lags behind the wave point, so that 
the phase velocity appears to be less than co. 

Now consider once more the estuary AF to which fig. 3 refers. The verti- 
cal tide H is supposed to be known as function of t at the inlet F and Q is 
known at the end A. Moreover we suppose that H and Q in the whole estuary 
at the instant tg are known as initial conditions or may be assumed as such. 

We divide the estuary into a number of sections of equal time of propaga- 
tion: 


TT) = VBM 


Let ABCDEF in fig. 3 represent such a division. Moreover the interval At 
of the time division is put equal to Tp. Then the sides of the lozenges are 
subcharacteristics. 

The computation by characteristics now proceeds as follows: From the 
states of motion associated with the tx-points 11 and 13 we deduce the state 
of motion associated with 22. Likewise we proceed from 13 and 15 to 24. 
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From 15 we proceed to 26 and define H there with the aid of the given value 
of Q. Then we go from 22 to 31 and use the boundary condition in F. Next we 
proceed from 22 and 24 to 33 and from 24 and 26 to 35. Then we go to 42 and 
soon. The order in which the points are treated might be chosen somewhat 
differently. 

If there were no or a negligible resistance the procedures represented in 
fig. 8 might be used. As the resistance in tidal motions is quite appreciable 
as a rule, we must modify the constructions as will be discussed below. Let 
P’ in fig. 10a represent the state of motion associated with the tx-point 43 in 
fig. 3 and let likewise Q' correspond to 41. The state of motion associated 
with 52 and represented by the HQ-point S', must be determined by construct- 
ing the contra-subcharacteristics associated with 43-52 and 41-52. 

Consider the wave point of which 43-52 represents the itinerary. It en- 
counters series of states of motion beginning with that represented by P’ and 
ending with that represented by S'. The factor |Q|Q on its journey varies 
from |Qpl Q, to |Q,|Qg- Let |Q7Q denote the average |Q |Q during this 
journey. Then we deduce from (509): 


(511) F, -Fp= - TOTO lt; - to) = 


Here Woz; is the resistance of the section DE. For co(t; - tg) is equal to the 
length 193 of that section and Wo, = wh 
Substituting from (503) in 6th yiel 


(512) a) (Hs - Hy + + - =0. 


This means that the point S' with the coordinates H, and Q, must lie on a line 

k sloping at Zo to 1 and through a point D which lies at the distance Wa, (Q7Q 

below P’. This line can be constructed when we make an estimation of |Q/@ 
In a similar manner we have 


(512) b) (Hg - Hg - Wy2!QIQ) - Z,(Q,- Qq) = 0 


and we can construct a line m in connection with the wave point of which 41-52 
is the itinerary. The point S' is then found at the intersection of k and m. 

The contra-subcharacteristics are not the lines k and m, but the lines P'S’ 
and Q'S’. 

In fig. 10a we have also represented the construction for the state of motion 
point R' associated with 67 (fig. 3). Here the datum Hp is introduced. 

The estimation of |Q|Q can be checked after the construction of S', and if 
necessary the construction is repeated with a new estimate. It appears, how- 
ever, that these estimations often require a relatively great amount of trial 
and error labour. In view of its practical importance we go somewhat deeper 
into the question: 

We approximate the variation of Q along the subcharacteristic 13-22 by 


Q = + - QM) 


for a moment adopting the place D as origin for x. Then, if Qp and Q, have 
equal signs, 


(513) 2,49 9, +19,2) 
tas 3% * 5% 


This is a function of Q, which is known and Q, which is unknown and ought to 
be estimated. Unfortunately § is very sensitive to errors in the estimation of 
Qs and this is the root of the difficulties expounded above. 

In order to substantially eliminate this drawback, we rewrite (513) as 
follows: 
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Qp + Qn (Qs - Q,) where = [2 a, 


The estimation of Q, is now only used to determine the factor Q,, which is 
relatively little sensitive. Then (512) may be reduced to 


(514) (H, - Hy +Z(Q,-Q)=0 » 
where = Wo3|Q, and Z = Zo + W239, 


Since Q,, is determined by estimation, (514) constitutes a linear equation in 
Q,. It is constructively interpreted in the way described already; the distance 
P'D becomes H,, and the slope of k becomes Z to 1. 

When Q) and have different signs, we may introduce the approximating 
formulae 


2 ‘ 1 


where we make use of an estimation for Q., and further proceed as described 
above. The justification of this procedure would demand a disproportionate 
space and hence is omitted here. 

A slightly different method with fixed subcharacteristics was applied by 
Lamoen to the Panama canal(25). The canal was divided into six sections. 
The results of Lamoen's computation are represented in fig. 6. 


5,3. Exact computation by variable velocities of propagation. 


We apply the characteristic transformation (cf (28) Ch 2) to the exact equa- 
tions (207) and (208). This means that we should add Z + times (207) to (208), 
where Z + is a factor to be chosen in such a way that we arrive at equations 
of the form 


(515) 1 08 , zt 
& ox c= 


Elaboration yields 


b - v2 
(516) a) =v V1 + =v-cg + 

b - v2 
(517) a) Z” = Zo and b) Z = -Zp Yi + = 


Here cp stands for 1/V¥bm and vi for ¥m/b as in (505); however, now m and b 
and hence co and Zp are not constant but depending on the state of motion, i.e. 
on Hand Q. Furthermore v denotes the velocity of flow v = Q/ and Vo is the 


critical velocity of flow defined by 
(518) Ve = VEATB, = coVb/b, , 
to which we return further below. 


It follows from (515) that in a wave point moving with the velocity dx /dt = 
c*, the states of motion satisfy the differential equation 


(519-a) + w{Q{Q dx = Z~ dQ [if dx=ctdt) 
Likewise 
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(519-b) dH + w |Q{Q dx = Z* dQ [if dx=c™ dt) 


holds good jn a wave point moving with the velocity c°. We might bring (519- 
a) or b) in a form like (509) by introducing the variables F and G, defined as 
functions of H and Q by 


OH OH 


However, the solution of these equations, in which Z*+ and Z~ are functions of 
H and Q and moreover of x, cannot be generally formulated. It moreover 
serves no practical purpose since we can as well operate directly with (519-a) 
and b). 

The velocity of propagation c* or c~ can be either positive or negative de- 
pending on the velocity of flow. The critical value of the velocity of flow is 
Vc defined by (518). When -v,<v<v_,, the flow is subcritical (flowing water) 
and c*+>0 and c~<0, When vsv¢ or v< -v,, the flow is supercritical (running 
or shooting water) and c+>c~>0 or c”“<c*<0. The distinction is of great prac- 
tical consequence for the influence of boundary conditions on the flow. We 
sha} ot go further into this question and confine ourselves to subcritical flow 

Ch 3, sect 212). 

Tidal motions usually are largely subcritical, i.e. v is small compared to 

V~- We may then be justified in approximating by 


(520) a)c*= veg and b) Z* = 


and by determining cg, Zo and w as functions of H by neglecting v2/2g. 

We now consider again the estuary AF mentioned before. A division into 
sections of approximately equal time of propagation is established. Further- 
more suppose for the moment that the tidal curves are smooth continuous 
functions of time. We can then perform the computation according to a grid 
as shown by fig. 3 in which the lozenges keep within the subcharacteristic tri- 
angles. We shall not describe this systematically but only give an illustrative 
example: 

Suppose the computation has been completed as far as the instant tg so that 
we know the states of motion associated with the tx-points 21, 23, 32, 41, 43 
etc. Now we proceed to 52 by considering the two wave points meeting each 
other in the place E at the instant t ,(tx-point S in fig. 10b). Let t, be the in- 
stant at which the descending wave S oint passed by D(tx-point P in fig. 10b) 
and th the instant at which the ascending wave point left F(tx-point Q in fig. 
10b). 

The subcharacteristic representing the itinerary of the descending wave 
point goes through Pand S. The straight line PS which is a chord of this sub- 
characteristic, is constructed as follows: 

Let (Hp, Q, ) and (H,, Q,) be the states of motion associated with P and S. 
Then we put 

lo and H, +25, 

and derive cg = (Qg/A) + Co (H,) from them as an approximation for the aver- 
age velocity of the descending Grave point. This will be sufficiently accurate 
tae the sections of the estuary are not too great. Now we estimate Hp, 

H, and Qs, compute cf, and construct P by drawing PS through S by the 
pg: Cg to 1. Then the estimations of H, and Qp can be checked by interpolat- 
ing betueen the states of motion in D at the instants to and t4 which were sup- 
posed to be known. The check of H, and Q, follows later. 
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In a similar way we construct SQ and interpolate for Hg and . Then the 
HQ-points P' and Q' associated with P and Q are known and S' is constructed 
as indicated in fig. 10a. This provides the check for the estimation of H, and 
Qs. If necessary the whole set of constructions is repeated. 

The interpolations mentioned may be performed graphically by two auxilia- 
ry diagrams, an Ht- and a Qt-diagram. 

We conclude this section by a few remarks: 

The approximations and estimations introduced above do not form an es- 
sential feature of the method, They usually meet the requirements in the 
Dutch practice, where we have channels of some 5 or 10 m depth and where 
we can operate with sections of 5 to 10 km with times of propagation of some 
10 or 15 min. In different situations modifications in the performance should 
be considered. The method could be applied practically unmodified to the 
Panama sea level canal which was divided into 4 sections each with a time of 
propagation of about 20 min. The results are shown in fig. 6. 

A fixed grid like that of fig. 3 can be used profitably when the tidal motion 
shows a gradual trend. When there are acute bends in the tidal curves, e.g. in 
case of manipulations with locks etc., the fixed grid must be abandoned and 
the characteristics marking the propagation of the sharp details should be 
followed uninterruptedly. 

Instead of graphical constructions, a numerical procedure may be used. 
This may be preferable when computing machines are used. From the en- 
gineer's point of view the graphical constructions have the advantage of help- 
ing to visualize the procedure. 


5,4. Jumps and bores, 


The bore is a tidal phenomenon observed in rather shallow estuaries and 
rivers with a great tidal range, which occurs in particular when the estuary 
or river is funnel-shaped. The bore comes into being when higher parts of 
the rising tide front of the penetrating tidal wave tend to overtake lower parts. 
This is put into evidence by the tx-diagram where the concurrent subcharac- 
teristics of the wave points in the rising tide front are intersecting. 

The bore is a hydraulic jump which is not fixed in place but travelling up 
the river. Similar mobile jumps may arise from too abrupt operations with 
locks etc. 

Although in these mobile jumps the vertical accelerations are of essential 
influence, it is still possible to compute a tide with a bore or other jump as a 
long wave, provided we account for the jump as a discontinuity in the long wave 
localized in a definite moving point comparable to a wave point. Doing so is 
justified since the length of the bore is as a rule very small as compared to 
the length of the tidal wave. 

A rigorous treatment of the mobile jump by rather simple formulae can be 
given if the cross-section of the channel is rectangular (cf (28) Ch 12 sect 31). 
The formulae for an arbitrary cross-section with a storing width b different 
from the width bg of the conveying streambed, are becoming unworkably in- 
volved, The following approximate formulae, however, may be used if the 
height of the jump is not too great, say less than half the depth: 

Let hj be the lower level ahead and hg the higher level in the rear ofa 
jump, which travels up a river in the negative sense of x. Let Hj, Q), Vp 
C91, Zo, etc. be associated with hy, and Hg etc. with hg. Then the velocity c 
of the jump is approximated by 


(521) c= (cg, - v4) +3 (egg - vg) 
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and between the states of motion separated by the jump a relation approximat- 
ed by 


(522) Hy - H, + Zo1 +} - = 0 


exists. 

The jump moves faster than the concurrent wave points ahead, but slower 
than the concurrent wave points in the rear. So both kinds of wave points 
meet the jump,and merge into it. Considered by an observer moving with the 
jump the flow in the low level region appears as supercritical and that in the 
high level region as subcritical. 

The jump is progressive, i.e. moving toward the low level region, if c>0. 
Such are the jumps in tidal regions as a rule. When c<0, the jump would be 
regressive. We have to deal with the well-known stationary jump if c = 0. 


Fig 11. Constructions to account for a 
jump (bore) ina tide 


In order to demonstrate how a jump is accounted for in a tidal computa- 
tion by characteristics, let fig. 11 represent part of such a computation. Sup- 
pose the jump at the instant of passing by the place C(tx-point S) has been de- 
termined. This is represented by the discontinuities S,Sg in the Ht-diagram 
and the Qt-diagram. We set out to determine the jump as it passes by the 
place B. To that purpose we consider three wave points, all meeting the jump 
in B. They are: a concurrent wave point running ahead (itinerary RP), a con- 
current wave point coming in the rear (itinerary TP), and a wave point en- 
countering the jump (itinerary QP). 
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First we compute the velocity c of the jump according to (521), using esti- 
mations where necessary, and we construct the tx-point P by drawing the 
line i: through S by the slope c to 1. 

Next the state of motion (HQ-point P}) just in front of the jump when it 
reaches B, is determined by considering the concurrent wave point ahead and 
the wave point coming from A. Then, according to (522) we draw a line m 
through the HQ-point P; by the slope 1/2 (Zo, + Zgg) tol. This line together 
with the construction associated with the wave point following the jump, de- 
fines the HQ-point Pg representing the state of motion just in the rear of the 
jump when it reaches B. 


CONCLUSIONS 


6,1. Comparative appreciation of computation methods 


As the rather great variety of the computation methods presented in this 
paper may appear somewhat bewildering, we shall now endeavour to give a 
comparitive appreciation. 

It is clear that, if there were a simply formulated exact solution of the 
mathematical equations of the tidal motion, there would be no need for any 
other solution. The tidal problems however are so intricate that such a solu- 
tion is not possible and so we must accept the existence of various methods 
approaching the problem in different ways. 

From the point of view of applicability we have to distinguish between ap- 
proximate and exact methods, By an exact method we understand a method 
by which the solution of the mathematical problems can be determined to any 
desired degree of accuracy. This means that the accuracy which can be ob- 
tained practically is limited only by the accuracy of the observational data 
from which the computation starts. 

The improvement of the accuracy generally goes at the cost of more 
labour. It depends on the purpose pursued how far one should go. For an ex- 
plorative investigation approximate computations generally will do. A more 
detailed investigation requires greater accuracy and then it depends for a good 
deal on the nature of the problem which method is the most appropriate. 

For the execution of the computation it makes a great difference whether 
trial and error procedures have been accepted in the method, or the computa- 
tion proceeds straight forward. Whereas a straight forward procedure can be 
entrusted to a relatively unexperienced computer, a trial and error procedure 
can only be efficient in the hands of a computer with great experience. 

We shall now first discuss the methods more in detail and then conclude by 
considering the use of computing machines. 


Harmonic methods. The single-harmonic method in many cases repro- 
duces very satisfactorily the fundamental of a periodic tide. The linearization 
of the resistance which has to be based on an estimation of the discharges, can 
be improved by successive approximations. An experienced computer often 
can make a fair estimate at once. Further the method is a straight forward 
procedure demanding relatively little labour. It is very appropriate for rapid 
exploration. 

The sine approximation of a tide is not always acceptable. The tidal flow 
in particular may deviate appreciably from the sine trend. Then one may 
adopt a procedure of successive approximations, such that every approxima- 
tion is extended to one more higher harmonic than the preceding approxima- 
tion. The first steps in this procedure generally improve the accuracy, but 
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the further steps give lesser improvements although the necessary labour 
increases. It is even dubious if the whole procedure is convergent at all. 
Practically it is usually not economic to go beyond the second harmonic. This 
double-harmonic method is still approximate. 

The application of a harmonic method to a complicated network of channels 
necessitates a special analysis to deal with the connections between the chan- 
nels, which may demand an appreciable amount of labour in excess to the 
computations of the separate channels. The network analysis is a straight 
forward procedure (cf Ch 4, sect 3). 

The harmonic methods are particularly well suited to deal with periodic 
tides. This restriction is not necessary in principle. For firstly we may 
treat nonperiodic functions by Fourier integrals, but the practical difficulties 
then encountered are in fact prohibitive. Secondly it has been endeavored to 
treat nonperiodic motions, in particular storm tides by suitable approximate 
functions which are readily integrable. Lorentz(12) tried to work on the as- 
sumption of periodically occurring storm surges, Mazure(17) put forward the 
approximation by an exponentially exploding sine function, and we might also 
try sums of real exponential functions. All these artifices can not lead but to 
approximate methods because an improvement of the accuracy along these 
lines requires excessively much labour. For accurate computations of non- 
periodic motions we may therefore exclude the harmonic methods. 


Direct methods. Among the methods of quad-differences, power series and 
iteration, the latter is the most refined. It lends itself very well for the nu- 
merical analysis of observed tidal motions, which forms the indispensible 
check of the schematization (cf 2,3). Therefore the iterative method is in 
particular efficient in estuaries which are hard to schematize. 

Prediction by one of the above methods requires the simultaneous solution 
of a number of nonlinear equations. In systems of a few sections this can still 
be done in a straight forward manner. In greater systems one must accept a 
trial and error procedure. In networks of simple structure this is still feasi- 
ble, but in complicated networks the trial and error labour rapidly becomes 
prohibitive. 

The method of cross-differences is substantially a straight forward pro- 
cedure. In Holsters’ original form it is an approximate method, reproducing 
more details than a single-harmonic method but also requiring more labour. 
For explorative computations it may be appropriate. As soon as the method 
has to be refined, the labour involved increases rather rapidly. 


Characteristic methods. The most profitable simplification in this type of 
method, viz the neglect of the resistance, is seldom admissable in tidal prob- 
lems. Computing the resistance demands rather much labour. Therefore the 
characteristic methods lend themselves not so much for explorative as for 
detailed investigations. 

In comparison to the iterative method, an exact computation along charac- 
teristics usually demands more labour when the tidal motion considered is an 
observed tide or not too much deviating from such tides, and the tidal system 
is not a too complicated network. One of the advantages of the characteristic 
method however, lies in the straight forward procedure which enabies us to 
predict tidal motions in complicated networks with an amount of labour rough- 
ly proportional to the extent of the system. 

Characteristic methods are particularly well suited to deal with waves of 
finite extent such as produced by locks, dam failures, etc. 
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Computing machines. The great amount of computing labour necessary to 
obtain accurate results on a tidal problem asks for considering the use of a 
computing machine, in particular of a rapid electronic computer. The ef- 
ficiency of such a machine depends largely on the possibility of arranging the 
computation in the form of a programme according to which the machine can 
work on uninterruptedly. Tidal computations however lend themselves to this 
only partly. 

According to Dutch experience a considerable part of the labour of a tidal 
computation has to be devoted to finding the most appropriate schematization 
of the tidal system. Here much depends on the judgment of the computer and 
for that reason a programme can not be well given. When a schooled comput- 
ing team has been formed to deal with this preliminary work, it can also deal 
with the proper predicting computations. Hence the use of a programme com- 
puting machine has only come into consideration in the Netherlands since very 
recently the task of the tidal hydraulicians was increased considerably. 

The use of a computing machine can hardly be expected to be economic 
for approximate computations. To set up a programme, a straight forward 
procedure is moreover requisite. Hence the methods of characteristics and 
cross-differences are the most promising in this way of approach. 


6,2. Comparative discussion of computations and model research . 


Model research and computation methods both pretend to yield solutions of 
tidal problems. A comparison therefore should not be omitted. 

There are mainly two types of models: hydraulic models and electric 
analogues, which we describe briefly: 


Hydraulic model. Geometrically true reduced scale models of extensive 
tidal systems should be made very large in order to observe certain limits to 
the scale reduction, firstly because the Reynolds' number in the model should 
be sufficiently great and secondly in view of the accuracy of measuring the 
vertical tide. As both conditions apply substantially to the vertical scale, 
models have been adopted in which the geometrical similarity is no longer 
strictly observed (distorted models). Such models are necessarily more or 
less schematized in a degree depending on the rate of distortion, and a true 
representation of the local flow patterns is no longer assured. Also the dis- 
tribution of velocity in a cross-section is affected, but it can be demonstrated 
that the total discharges are truly represented provided the resistance is in- 
creased in an appropriate manner (exaggerated roughness). In fact this ought 
to be checked section by section. 

A hydraulic model visualizes the water motion very clearly and directly 
which may be a great help for the engineer. As a disadvantage of hydraulic 
models may be noted that it is often difficult (partly because of the exaggerat- 
ed roughness) to measure the variable discharges. 


Electric model. On the basis of the analogy of electric and hydraulic 
systems(9°) it is in principle possible to make an electric analogue of a tidal 
system(26), Between the principle and the practical execution there is a long 
way. The model is made section-wise so that each model section represents 
a corresponding channel section in such a way that the total storage, inertance 
and resistance are truly represented. The internal mechanism in the model 
section offers no analogy with the channel section. The electric model there- 
fore only represents truly the total discharges, which means that it is a 
schematic model, The electric currents analogous to these discharges can 
in principle be recorded just like the potentials representing the heads. 
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An electric model does not visualize the tidal motion as a hydraulic model 
does. It offers great possibilities however to produce very rapidly visual 
records of a great variety of tidal diagrams. 

The choice between computations and models may be governed by the fol- 
lowing considerations: 


1. Typical advantages. Models offer rapid visualization possibilities. On 
the other hand computations enable us to penetrate more deeply into the physi- 
cal mechanism of the motion, thus improving the insight. This is in particu- 
lar true for computations with graphs and slide rules, and in a much lesser 
degree for computations with programme computing machines. 

It depends partly on the nature of a tidal problem, whether a model or com- 
putations will be most appropriate, and often both means may profitably be 
applied in cooperation, 


2. Availability. Not always the facilities for both computing and model 
research are at hand. Model research requires a well equipped laboratory 
and computations require trained computers. 


3. Accuracy. The accuracy of both models and computations is limited by 
the errors of many preliminary data, such as the Chézy coefficient. Also 
computations as well as distorted hydraulic and electric models are to a cer- 
tain extent dependent on schematization. A generally valid comparison as to 
the degree of the accuracy between the three seems hardly possible as the 
means by which it can be tried to reduce the deviations from nature are not 
identical. 

In hydraulic models both the effects of schematization and the measuring 
errors can be decreased by the use of larger models. Here the economy be- 
comes a very important argument (v. below). 

In an electric model the accuracy depends on the schematization and partly 
on the veracity of the representation of hydraulic properties (such as the quad- 
ratic resistance) by electric elements. Whereas the former can be improved 
rather easily by increasing the number of sections by which the schematiza- 
tion is set up, the latter can only be improved as far as the electric technique 
admits. 

The accuracy of a computation can be improved by refining the schematiza- 
tion and, if we consider exact methods, by proceeding to a further approxima- 
tion. 


4. Economy. There is a rather great difference between models and com- 
putations in the ratio of variable to fixed costs. The building of a model is a 
rather expensive affair. Once the model built however, it is relatively easy 
to deal with a great variety of problems in the same tidal system. The prep- 
aration for a computation (schematization, checking, etc.) likewise takes much 
labour, although relatively less than in case of a model, whereas the investi- 
gation of a number of problems takes relatively more labour than in a model. 

In the Netherlands both types of models as well as several of the methods 
of computation treated in the preceding chapters, play a useful part in the 
solution of the manifold problems with which the tides confront us. 
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